% Paper: E- Monetary Theory  (2nd edition)
% Author: Duong Ngotran
% Graph: Crisis - Shock on kappa - Taylor Rule

%----------------------------------------------------------------------------
% Some Parameters
%----------------------------------------------------------------------------
T = 301;    %Number of periods 
var_foc = 25;      % Number of variables in one period
pos_value = 3;     % Position of var_value in A
truc = (-1):(T-2);         % Horizontal axis when drawing with level
spline_truc = (-1):(1/4):(T-2);     %Horizontal axis when drawing with spline 
size_struc = 1201;  % size of spline truc
deltab = 0.5;

%----------------------------------------------------------------------------
% Mapping name of variables to its position
%----------------------------------------------------------------------------
mapp = containers.Map;
mapp('gamma')=1; mapp('cb')=2; mapp('Rf')=3; mapp('Rm')=4; mapp('ql')=5;
mapp('mur')=6; mapp('muc')=7; mapp('ip')=8; mapp('pie')=9; mapp('n')=10;
mapp('m')=11; mapp('bh')=12; mapp('k')=13; mapp('tau')=14; mapp('s')=15;
mapp('ch')=16; mapp('etaz')=17; mapp('lama')=18; mapp('lamb')=19; mapp('etab')=20;
mapp('i')=21; mapp('pm')=22; mapp('y')=23; mapp('l')=24; mapp('u')=25;

%----------------------------------------------------------------------------
%Read data file into matrix A (time X var_position X var_value) 
%----------------------------------------------------------------------------
A = zeros(T*var_foc, pos_value);
A = dlmread('EMT_ir_shock.txt');

%----------------------------------------------------------------------------
%Create matrix AT: row (time)X col (var_position)
%----------------------------------------------------------------------------
AT = zeros(T, var_foc);
for t=1:T
    for i=1:var_foc
        AT(t,i)= A((t-1)*var_foc+i, pos_value);
    end;
end;

%----------------------------------------------------------------------------
%Create matrix AP: Percentage deviation from steady state
%AP contains 2 columns: time X percentage_deviation
%----------------------------------------------------------------------------
% Vector stores the steady state "ss"
ss = zeros(1, var_foc);
for i=1:var_foc
    ss(1,i) = AT(T,i);
end;

% Matrix AP
AP = zeros(T, var_foc);
for t=1:T
    for i=1:var_foc
        AP(t,i)= (AT(t,i)-ss(i))/ss(i)*100;
    end;
end;
        
%----------------------------------------------------------------------------
%Create matrix AP_S: Percentage deviation from steady state (spline)
%AP contains 2 columns: time X percentage_deviation
%----------------------------------------------------------------------------       
AP_S = zeros(var_foc, size_struc);
for i=1:var_foc
    s_var = spline(truc, AP(:,i), spline_truc);
    AP_S(i,:)= spline(truc, AP(:,i), spline_truc);
end;

%----------------------------------------------------------------------------
%Graph 1: Federal Funds rate (level0, real FFR, real lending rate (anually)
%---------------------------------------------------------------------------- 
rf_t = zeros(1,T);   % Federal Funds Rate Level
rrf_t = zeros(1,T);  % Real Federal Funds Rate
rh_t = zeros(1,T);   % Real Borrowing rate 
for t=1:T
    rf_t(t) = (AT(t,mapp('Rf'))-1)*400;
end;
for t=1:(T-1)
    rrf_t(t) = (AT(t,mapp('Rf'))/AT(t+1,mapp('pie'))-1)*400;
    rh_t(t) = ((deltab+ deltab*AT(t+1,mapp('ql')))/AT(t,mapp('ql'))/AT(t,mapp('pie'))-1)*400;
end;
rrf_t(T) = rrf_t(T-1);
rh_t(T) = rh_t(T-1);
rrf_t(1) = rrf_t(T);   % Before the unexpected shock
rh_t(1) = rh_t(T);

%----------------------------------------------------------------------------
%Graph 2: C_bankers vs C_household vs Aggregate C
%----------------------------------------------------------------------------
agg_c = zeros(T, 1);   %Aggregate consumption
percent_agg_c = zeros(T, 1); % Percentage deviation of aggregate consumption from steady state

for t=1:T
    agg_c(t) = AT(t,mapp('cb'))+ AT(t,mapp('ch'));
end;
agg_c_ss = agg_c(T); % Steady state of aggregate consumption
for t=1:T
    percent_agg_c(t) = (agg_c(t)-agg_c(T))/agg_c(T)*100;
end;
s_agg_c = spline(truc, percent_agg_c, spline_truc);

%----------------------------------------------------------------------------
%Graph 4: Nominal growth M(t)/M(t-1), N(t)/N(t-1), B(t)/B(t-1)
%----------------------------------------------------------------------------
growth_M = zeros(1,T);
growth_B = zeros(1,T);
growth_N = zeros(1,T);
for t=2:T
    growth_M(t) = (AT(t,mapp('m'))/AT(t-1,mapp('m'))*AT(t,mapp('pie'))-1)*100;
    growth_B(t) = (AT(t,mapp('bh'))/AT(t-1,mapp('bh'))*AT(t,mapp('pie'))-1)*100;
    growth_N(t) = (AT(t,mapp('n'))/AT(t-1,mapp('n'))*AT(t,mapp('pie'))-1)*100;
end;

s_growth_M = spline(truc, growth_M, spline_truc);
s_growth_B = spline(truc, growth_B, spline_truc);
s_growth_N = spline(truc, growth_N, spline_truc);


%----------------------------------------------------------------------------
%List of Figure 
%----------------------------------------------------------------------------
set(groot, 'DefaultLegendInterpreter', 'latex')

figure('Name', 'Federal Funds Rate and Real Borrowing Rate');
plot(truc, rf_t, truc, rh_t);
leg = legend('$R^f$ - Taylor', '$r^h$ - real borrowing rate', 'location', 'best'); 

figure('Name', 'Real balance of ZMDs and Reserves');
plot(spline_truc, AP_S(mapp('m'),:), spline_truc, AP_S(mapp('n'),:));
leg = legend('$m$ - ZMD', '$n$ - Reserves', 'location', 'best');

figure('Name', 'Consumption of bankers and households');
plot(spline_truc, AP_S(mapp('cb'),:), spline_truc, AP_S(mapp('ch'),:));
leg = legend('$c$ - Banker consumption', '$\tilde{c}$ - Household consumption', 'location', 'best');

figure('Name', 'Investment and Capital');
plot(spline_truc, AP_S(mapp('i'),:), spline_truc, AP_S(mapp('k'),:));
leg = legend('$i$', '$k$', 'location', 'best');

figure('Name', 'Inflation');
plot(spline_truc, AP_S(mapp('pie'),:));
leg = legend('$\pi$', 'location', 'best');

figure('Name', 'Output');
plot(spline_truc, AP_S(mapp('y'),:));
leg = legend('$y$', 'location', 'best');
        
        
        
        
        
        
        